A variational eigenvalue solver on a quantum processor 
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Quantum computers promise to efficiently solve important problems that are intractable on a con- 
ventional computer. For quantum systems, where the dimension of the problem space grows expo- 
nentially, finding the eigenvalues of certain operators is one such intractable problem and remains a 
fundamental challenge. The quantum phase estimation algorithm can efficiently find the eigenvalue 
of a given eigenvector but requires fully coherent evolution. We present an alternative approach 
that greatly reduces the requirements for coherent evolution and we combine this method with a 
new approach to state preparation based on ansatze and classical optimization. We have imple- 
mented the algorithm by combining a small-scale photonic quantum processor with a conventional 
computer. We experimentally demonstrate the feasibility of this approach with an example from 
quantum chemistry — calculating the ground state molecular energy for He-H+, to within chem- 
ical accuracy. The proposed approach, by drastically reducing the coherence time requirements, 
enhances the potential of the quantum resources available today and in the near future. 
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In chemistry, the properties of atoms and molecules can 
be determined by solving the Schrodinger equation. How- 
ever, because the dimension of the problem grows expo- 
nentially with the size of the physical system under con- 
sideration, exact treatment of these problems remains 
classically infeasible for compounds with more than 2-3 
atoms [1^. Many approximate methods [1 have been de- 
veloped to treat these systems, but efficient exact meth- 
ods for large chemical problems remain out of reach for 
classical computers. Beyond chemistry, the solution of 
large eigenvalue problems [3 would have applications 
ranging from determining the results of internet search 
engines [4 to designing new materials and drugs [5j. 

Recent developments in the field of quantum compu- 
tation offer a way forward for efficient solutions of many 
instances of large eigenvalue problems which are classi- 
cally intractable [6HT2]. Quantum approaches to finding 
eigenvalues have previously relied on the quantum phase 
estimation (QPE) algorithm. The QPE algorithm offers 
an exponential speedup over classical methods and re- 
quires a number of quantum operations 0{l/p) to obtain 
an estimate with precision p pISHIS] . In the standard 
formulation of QPE, one assumes the eigenvector |?/;) of 
a Hermitian operator H is given as input and the prob- 
lem is to determine the corresponding eigenvalue A. The 
time the quantum computer must remain coherent is de- 
termined by the necessity of 0{l/p) successive applica- 
tions of e""^^^, each of which can require on the order of 
millions or billions of quantum gates for practical appli- 
cations [171 [E] 5 as compared to the tens to hundreds of 
gates achievable in the short term. 

Here we introduce and experimentally demonstrate an 
alternative to QPE that significantly reduces the require- 



ments for coherent evolution. We have developed a re- 
configurable quantum processing unit (QPU), which effi- 
ciently calculates the expectation value of a Hamiltonian 
(H), providing an exponential speedup over conventional 
methods. The QPU is combined with an optimization al- 
gorithm run on a classical processing unit (CPU), which 
variationally computes the eigenvalues and eigenvectors 
of H. By using a variational algorithm, this approach 
reduces the requirement for coherent evolution of the 
quantum state, making more efficient use of quantum 
resources, and may offer an alternative route to practical 
quantum-enhanced computation. 

Algorithm 1: Quantum expectation estimation 
This algorithm computes the expectation value of a given 
Hamiltonian Ti for an input state |?/^). Any Hamiltonian 
may be written as 
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for real h where Roman indices identify the subspace on 
which the operator acts, and Greek indices identify the 
Pauli operator, e.g. a = x. By exploiting the linearity of 
quantum observables, it follows that 



in) = ^hl{al,) + ^h%{ai,a^^) 
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We consider Hamiltonians that can be written as a num- 
ber of terms which is polynomial in the size of the sys- 
tem. This class of Hamiltonians encompasses a wide 
range of physical systems, including the electronic struc- 
ture Hamiltonian of quantum chemistry, the quantum 
Ising Model, the Heisenberg Model [20, 21j, matrices that 
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FIG. 1. Architecture of the quantum- variational eigensolver. Algorithm 1: Quantum states that have been previously 
prepared, are fed into the quantum modules which compute (Hi), where Tii is any given term in the sum defining Ti. The 
results are passed to the CPU which computes (Ti). Algorithm 2: The classical minimization algorithm, run on the CPU, 
takes {H) and determines the new state parameters, which are then fed back to the QPU. 



are well approximated as a sum of n-fold tensor prod- 
ucts mill], and more generally any /c— sparse Hamilto- 
nian without evident tensor product structure (see Ap- 
pendix for details). Thus the evaluation of (H) reduces 
to the sum of a polynomial number of expectation val- 
ues of simple Pauli operators for a quantum state |?/^), 
multiplied by some real constants. A quantum device 
can efficiently evaluate the expectation value of a tensor 
product of an arbitrary number of simple Pauli opera- 
tors [5 , therefore with an n-qubit state we can efficiently 
evaluate the expectation value of this 2"^ x 2^ Hamilto- 
nian. 

One might attempt this using a classical computer 
by separately optimizing all reduced states correspond- 
ing to the desired terms in the Hamiltonian, but this 
would suffer from the A^-representability problem, which 
is known to be intractable for both classical and quantum 
computers (it is in the quantum complexity class QMA- 
Hard [24] ) . The power of our approach derives from the 
fact that quantum hardware can store a global quantum 
state with exponentially fewer resources than required by 
classical hardware, and as a result the N-representability 
problem does not arise. 

As the expectation value of a tensor product of an arbi- 
trary number of Pauli operators can be measured in con- 
stant time and the spectrum of each of these operators is 
bounded, to obtain an estimate with precision p, our ap- 
proach incurs a cost of 0{\h\'^/p'^) repetitions. Thus the 
total cost of computing the expectation value of a state 
IV^) is given by 0{\hmax\'^M/p'^)^ where M is the number 
of terms in the decomposition of the Hamiltonian. The 
advantage of this approach is that the coherence time to 
make a single measurement after preparing the state is 



0(1). In essence, we dramatically reduce the coherence 
time requirement while maintaining an exponential ad- 
vantage over the classical case, by adding a polynomial 
number of repetitions with respect to QPE. 
Algorithm 2: Quantum variational eigensolver 
The procedure outlined above replaces the long coher- 
ent evolution required by QPE by many short coherent 
evolutions. In both QPE and Algorithm 1 we require a 
good approximation to the ground state wavefunction to 
compute the ground state eigenvalue and we now consider 
this problem. Previous approaches have proposed to pre- 
pare ground states by adiabatic evolution [15 , or by the 
quantum metropolis algorithm [25 . Unfortunately both 
of these require long coherent evolution. Algorithm 2 is a 
variational method to prepare the eigenstate and, by ex- 
ploiting Algorithm 1, requires short coherent evolution. 
Algorithm 1 and 2 and their relationship are shown in 
Fig. [1] and detailed in the Appendix. 

It is well known that the eigenvalue problem for an 
observable represented by an operator H can be re- 
stated as a variational problem on the Rayleigh-Ritz quo- 
tient [26l[27], such that the eigenvector |?/;) corresponding 
to the lowest eigenvalue is the Itp) that minimizes 
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By varying the experimental parameters in the prepa- 
ration of IV^) and computing the Rayleigh-Ritz quotient 
using Algorithm 1 as a subroutine in a classical mini- 
mization, one may prepare unknown eigenvectors. At 
the termination of the algorithm, a simple prescription 
for the reconstruction of the eigenvector is stored in the 
final set of experimental parameters that define l?/^). 
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FIG. 2. Experimental implementation of our scheme, (a) Quantum state preparation and measurement of the expectation 
values {^l^\(Ji cTjlV^) are performed using a quantum photonic chip. Photon pairs, generated using spontaneous parametric 
down-conversion are injected into the waveguides encoding the |00) state. The state |'0) is prepared using thermal phase shifters 
01-8 (orange rectangles) and one CNOT gate and measured using photon detectors. Coincidence count rates from the detectors 
Di_4 are passed to the CPU running the optimization algorithm. This computes the set of parameters for the next state and 
writes them to the quantum device, (b) A photograph of the QPU. 



If a quantum state is characterized by an exponentially 
large number of parameters, it cannot be prepared with 
a polynomial number of operations. The set of efficiently 
preparable states are therefore characterized by polyno- 
mially many parameters, and we choose a particular set 
of ansatz states of this type. Under these conditions, a 
classical search algorithm on the experimental parame- 
ters which define l?/^), needs only explore a polynomial 
number of dimensions — a requirement for the search to 
be efficient. 

One example of a quantum state parametrized by a 
polynomial number of parameters is the unitary coupled 
cluster ansatz j^ 
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where T is the cluster operator (defined in the Appendix) 
and |^)^g^ is some reference state, normally taken to be 
the Hartree-Fock ground state. There is currently no 
known efficient classical algorithm based on these ansatz 
states. However, non-unitary coupled cluster ansatz is 
sometimes referred to as the "gold standard of quantum 
chemistry" as it is the standard of accuracy to which 
other methods in quantum chemistry are often compared. 



The unitary version of this ansatz is thought to yield 
superior results to even this "gold standard" [4 . Details 
of efficient construction of the unitary coupled cluster 
state using a quantum device are given in the Appendix 
(see also Ref. [29 ). 
Prototype demonstration 

We have implemented the QPU using integrated quan- 
tum photonics technology [30 . Our device, shown 
schematically in Fig. [2] is a reconfigurable waveguide chip 
that implements several single qubit rotations and one 
two-qubit entangling gate and can prepare an arbitrary 
two-qubit pure state. This device operates across the 
full space of possible configurations with mean statisti- 
cal fidelity F > 99% [31 . The state is prepared, and 
measured in the Pauli basis, by setting 8 voltage driven 
phase shifters and counting photon detection events with 
silicon single photon detectors. 

The ability to prepare an arbitrary two-qubit separable 
or entangled state enables us to investigate 4x4 Hamil- 
tonians. For the experimental demonstration of our al- 
gorithm we choose a problem from quantum chemistry, 
namely determining the bond dissociation curve of the 
molecule He-H+ in a minimal basis. The full configura- 
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FIG. 3. Finding the ground state of He-H^ for a specific 
molecular separation, i? = 90 pm. (a) Experimentally com- 
puted energy {V) (colored dots) as a function of the opti- 
mization step j. The color represents the tangle (degree of 
entanglement) of the physical state, estimated directly from 
the state parameters {0;^}. The red lines indicate the en- 
ergy levels of 1-L{R). The optimization algorithm clearly con- 
verges to the ground state of the molecule, which has small 
but non zero tangle. The crosses show the energy calculated 
at each experimental step, assuming an ideal quantum device. 
(b) Overlap | (V^-^IV^*^) | between the experimentally computed 
state l'^-^) at each the optimization step j and the theoretical 
ground state of 'H, IV^*^). Further details are provided in the 
Appendix. 



tion interaction Hamiltonian for this system has dimen- 
sion 4, and can be written compactly as 



u{R) = Y,Ki.R)< 



E h%{R)<^>- 



(5) 
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The coefficients hl^{R) and h^^JR) were determined us- 
ing the PSI3 computational package [3] and tabulated in 
the Appendix. 

In order to compute the bond dissociation of the 
molecule, we use Algorithm 2 to compute its ground state 
for a range of values of the nuclear separation R. In 
Fig. |3] we report a representative optimization run for a 
particular nuclear separation, demonstrating the conver- 
gence of our algorithm to the ground state oil-L{R) in the 
presence of experimental noise. Fig. IsFa) demonstrates 
the convergence of the average energy, while Fig. IsFb) 
demonstrates the convergence of the overlap | {ip^\iJP) \ 
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FIG. 4. Bond dissociation curve of the He-H^ molecule. This 
curve is obtained by repeated computation of the ground state 
energy (as shown in Fig. [sl for several 1-L{R). The magnified 
plot shows that after correction for the measured systematic 
error, the data overlap with the theoretical energy curve and 
importantly we can resolve the molecular separation of min- 
imal energy. Error bars show the standard deviation of the 
computed energy. 



of the current state \il)^) with the target state IV^^). The 
color of each entry in Fig. [sja) represents the tangle (ab- 
solute concurrence squared) of the state at that step of 
the algorithm. It is known that the volume of separable 
states is doubly-exponentially small with respect to the 
rest of state space [33 . Thus, the ability to traverse non- 
separable state space increases the number of paths by 
which the algorithm can converge and will be a require- 
ment for future large-scale implementations. Moreover, 
it is clear that the ability to produce entangled states is 
a necessity for the accurate description of general quan- 
tum systems where eigenstates may be non-separable, for 
example the ground state of the He-H+ Hamiltonian has 
small but not negligible tangle. 

Repeating this procedure for several values of i?, we 
obtain the bond dissociation curve which is reported in 
Fig. [4j This allows for the determination of the equilib- 
rium bond length of the molecule, which was found to be 
i?=92.3ib0.1 pm with a corresponding ground state elec- 
tronic energy of E= -2.865±0.008 MJ/mol. This energy 
has been corrected for experimental error using a method 
fully described in the Appendix. The corresponding theo- 
retical curve shows the numerically exact energy derived 
from a full configuration interaction calculation of the 



molecular system in the same basis. More than 96% of 
the experimental data are within chemical accuracy with 
respect to the theoretical values. At the conclusion of 
the optimization, we retain full knowledge of the exper- 
imental parameters, which can be used for efficient re- 
construction of the state |?/^) in the event that additional 
physical or chemical properties are required. 
Discussion 

Algorithm 1 uses relatively few quantum resources com- 
pared to QPE. Broadly speaking, QPE requires a large 
number of n-qubit quantum controlled operations to be 
performed in series — placing considerable demands on 
the number of components and coherence time — while 
the inherent parallelism of our scheme enables a small 
number of n-qubit gates to be exploited many times, 
drastically reducing these demands. Moreover, adding 
control to arbitrary unitary operations in practice is dif- 
ficult if not impossible for current quantum architectures 
(although a proposed scheme to add control to arbitrary 
unitary operations has recently been demonstrated [34 ). 
To give a numerical example, the QPE circuit for a 4 
X 4 Hamiltonian such as that demonstrated here would 
require at least 12 CNOT gates, while our method only 
requires one. 

In implementing Algorithm 2, the device prepares 
ansatz states that are defined by a polynomial set of 
parameters. This ansatz might be chosen based on 
knowledge of the physical system of interest (as for the 
unitary coupled cluster and typical quantum chemistry 
ansatze) thus determining the device design. However, 
our architecture allows for an alternative, and poten- 
tially more promising approach, where the device is first 
constructed based on the available resources and we de- 
fine the set of states that the device can prepare as the 
"device ansatz" . Due to the quantum nature of the de- 
vice, this ansatz can be very distinct from those used 
in traditional quantum chemistry. With this alterna- 
tive approach the physical implementation is then given 
by a known sequence of quantum operations with ad- 
justable parameters — determined at the construction of 
the device — with a maximum depth fixed by the coher- 
ence time of the physical qubits. This approach, while ap- 
proximate, provides a variationally optimal solution for 
the given quantum resources and may still be able to pro- 
vide qualitatively correct solutions, just as approximate 
methods do in traditional quantum chemistry (for exam- 
ple Hartree Fock). The unitary coupled cluster ansatz 
(Eq. [4| provides a concrete example where our approach 
provides an exponential advantage over known classical 
techniques. 

We have developed and experimentally implemented 
a new approach to solving the eigenvalue problem with 
quantum hardware. Algorithm 1 shares with QPE the 
need to prepare a good approximation to the ground 
state, but replaces a single long coherent evolution by 
a number of shorter coherent calculations proportional 



to the number of terms in the Hamiltonian. While the 
effect of errors on each of these calculations is the same 
as in QPE, the reliance on a number of separate calcula- 
tions makes the algorithm sensitive to variations in state 
preparation between the separate quantum calculations. 
This effect requires further investigation. In Algorithm 2, 
we experimentally implemented a ground state prepara- 
tion procedure through a direct variational algorithm on 
the control parameters of the quantum hardware. Larger 
calculations will require a choice of ansatz, for which 
there are two possibilities. One could experimentally im- 
plement chemically motivated ansatze such as the uni- 
tary coupled cluster method described in the Appendix. 
Alternatively one could pursue those ansatze that are 
most easy to implement experimentally — creating a new 
set of device ansatze states which would require classi- 
fication in terms of their overlap with chemical ground 
states. Such a classification would be a good way to 
determine the value of a given experimental advance — 
for ground state problems it is best to focus limited 
experimental resources on those efforts that will most 
enhance the overlap of preparable states with chemical 
ground states. In addition to the above issues, which we 
leave to future work, an interesting avenue of research is 
to ask whether the conceptual approach described here 
could be used to address other intractable problems with 
quantum-enhanced computation. Examples that can be 
mapped to the ground state problem, and where the n- 
represent ability problem does not occur, include search 
engine optimisation and image recognition. It should 
be noted that the approach presented here requires no 
control or auxiliary qubits, relying only on measurement 
techniques that are already well established. For exam- 
ple, in the two qubit case, these measurements are iden- 
tical to those performed in Bell inequality experiments. 
Quantum simulators with only a few tens of qubits are 
expected to outperform the capabilities of conventional 
computers, not including open questions regarding fault 
tolerance and errors/precision. Our scheme would allow 
such devices to be implemented using dramatically less 
resources than the current best known approach. 
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APPENDIX 

SUPPLEMETARY THEORY 

Quantum eigenvector preparation algorithm 

Below we detail the steps involved in implementing 
Algorithm 2. 

1. Design a quantum circuit, controlled by a set of 
experimental parameters {^i}, which can prepare 
a class of states. Using this device, prepare the 
initial state {ip^) and define the objective function 
im}) = {^m})\^\i'm})), which efficiently 
maps the set of experimental parameters to the ex- 
pectation value of the Hamiltonian and is computed 
in this work by Algorithm 1. n denotes the current 
iteration of the algorithm. 

2. Let n = 

3. Repeat until optimization is completed 

(a) Call Algorithm 1 with {Oi} as input: 

i. Using the QPU, compute (crj^,), {crl^cri), 
(al^aia^), ..., on {ip^) for all terms of H. 

ii. Classically sum on CPU the values from 
the QPU with their appropriate weights, 
h, to obtain /{{Of}) 

(b) Feed /({^f }) to the classical minimization 
algorithm (e.g. gradient descent or Nelder- 
Mead Simplex method), and allow it to deter- 
mine {<9[^+^}. 

Second Quantized Hamiltonian 

When taken with the Born-Oppenheimer approximation, 
the Hamiltonian of an electronic system can be generally 
written ll as 



^W = ^ hpq{R)alaq + ^/ip 



s{R)dlalaras 



(6) 



where a] and dj are the fermionic creation and annihila- 
tion operators that act on the single particle basis func- 
tions chosen to represent the electronic system and obey 
the canonical anti-commutation relations {dl^dj} = Sij 
and {di^dj} = {dj^d^A = 0. i? is a vector representing 
the positions of the Nuclei in the system, and is fixed for 
any given geometry. The constants hpq{R) and hpqrs{R) 
are evaluated using an initial Hartree-Fock calculation 
and relate the second quantized Hamiltonian to the first 
quantized Hamiltonian. They are calculated as 



dr Xpi^T 



h„„ = / *, *2 ^('•■)'x.('-2)-x,(n)x.(.-3) ,8) 
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where Xp(^) ^^^ single particle spin orbitals, Z^ is the 
nuclear charge, and r^ is the nuclear position. From the 
definition of the Hamiltonian, it is clear that the number 
of terms in the Hamiltonian is 0{N'^) in general, where N 
is the number of single particle basis functions used. The 
map from the Fermionic algebra of the second quantized 
Hamiltonian to the distinguishable spin algebra of qubits 
is given by the Jordan- Wigner transformation [2^, which 
for our purposes can be concisely written as 



«] 



) (J+ (8) cr 
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(9) 
(10) 



where (j+ and a- are the Pauli spin raising and lowering 
operators respectively. It is clear that this transforma- 
tion does not increase the number of terms present in the 
Hamiltonian, it merely changes their form and the spaces 
on which they act. Thus the requirement that the Hamil- 
tonian is a sum of polynomially many products of Pauli 
operators is satisfied. As a result, the expectation value 
of any second quantized chemistry Hamiltonian can be 
efficiently measured with our scheme. 

For the specific case of He-H+ in a minimal, STO- 
3G basis, it turns out that full configuration interaction 
(FCI) Hamiltonian has dimension four, thus a more com- 
pact representation is possible than in the general case. 
In this case, the FCI Hamiltonian can be written down 
for each geometry expanded in terms of the tensor prod- 
ucts of two Pauli operators. Thus the Hamiltonian is 
given explicitly by an FCI calculation in the PSI3 com- 
putational package [3 and can be written as 



niR) = Y^ humeri + J2 h%{R)<4 
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Unitary Coupled Cluster Theory 

One example of a state which is efficiently preparable 
on a quantum computer, but not so on a classical com- 
puter is the unitary coupled cluster expansion [4]. The 
unitary coupled cluster theory method is a variational 
ansatz which takes the form 



I*) 
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where |^)^g/' is some reference state, usually the Hartree 
Fock ground state, and T is the cluster operator for an 
N electron system defined by 



with 
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where repeated indices imply summation as in the main 
text, and higher order terms fohow logicahy. It is clear 
that by construction the operator (T — T^) is anti- 
hermitian, and exponentiation maps it to a unitary op- 
erator U = e^^~^ \ For any fixed excitation level /c, the 
reduced cluster operator is written as 



TW = ^ri 



(16) 



Unfortunately, in general no efficient implementation of 
this ansatz has yet been developed for a classical com- 
puter, even for low order cluster operators due to the 
non-truncation of the BCH series [4]. The reduced anti- 
hermitian cluster operator (T^^^ — T^^^''") is the sum of 
a polynomial number of terms in the number of one 
electron basis functions, namely it contains a number of 
terms 0{N^{M — N)^) where M is the number of sin- 
gle particle orbitals. By defining an effective Hermitian 
Hamiltonian H — i(T^^^ — T^^"^^) and performing the 
Jordan- Wigner transformation to reach a Hamiltonian 
that acts on the space of qubits, H, we are left with a 
Hamiltonian which is a sum of polynomially many prod- 
ucts of Pauli operators. The problem then reduces to 
the quantum simulation of this effective Hamiltonian, H, 
which can be done in polynomial time using the proce- 
dure outlined by Ortiz et al. |5]. This represents one 
example of a state which can be efficiently prepared on a 
quantum device, which cannot be efficiently prepared by 
any known means on a classical computer. 

Finding excited states 

Frequently, one may be interested in eigenvectors and 
eigenvalues related to excited states (interior eigenval- 
ues). Fortunately our scheme can be used with only mi- 
nor modification to find these excited states by repeating 
the procedure on Hx = (H — A)^. The folded spectrum 
method |7] allows a variational method to converge 
to the eigenvector closest to the shift parameter A. By 
scanning through a range of A values, one can recover 
the eigenvectors and eigenvalues of interest. Although 
this operation incurs a small polynomial overhead — the 
number of terms in the Hamiltonian is quadratic with 
respect to the original Hamiltonian — this extra cost is 
marginal compared to the cost of solving the problem 
classically. 

Application to /c— sparse Hamiltonians 

The method described in the main body of this work 
may be applied to general /c— sparse Hamiltonian ma- 
trices which are row-computable even when no efficient 
tensor decomposition is evident with only minor modi- 
ffcation. A Hamiltonian H is referred to as /c— sparse if 
there are at most k non-zero elements in each row and 
column of the matrix and row computable if there is an 



efficient algorthim for ffnding the locations and values of 
the non-zero matrix elements in each row of Ti. 

Let H be a 2"^ X 2^ /c— sparse row-computable Hamil- 
tonian. A result by Berry et al. [8 shows that H may be 
decomposed as H = ^T=i^j with m = 6/c^, Hj being 
a 1— sparse matrix and each Hj may be efficiently simu- 
lated (e~*^^^ may be acted on the qubits) by making only 
0(log* n) queries to the Hamiltonian H. Alternatively, a 
more recent result by Childs et al. [9 has found that it 
possible to use a star decomposition of the Hamiltonian 
such that m = 6k and each Hj is now a galaxy which can 
be efficiently simulated using 0(/c+log* N) queries to the 
Hamiltonian. Either of these schemes may be used to im- 
plement our algorithm efficiently for a general /c— sparse 
matrix, and the choice may be allowed to depend on the 
particular setup available. Following a prescription by 
Knill et al. [10 , the ability to simulate Hj is sufficient 
for efficient measurement of the expectation value (Hj). 
After determining these values, one may proceed as be- 
fore in the algorithm as outlined in the main text and 
use them to determine new parameters for the classical 
minimization. 



Classical optimization algorithm 

For the classical optimization step of our integrated 
processor we implemented the Nelder-Mead (NM) algo- 
rithm [11 , a simplex-based direct search (DS) method 
for unconstrained minimization of objective functions. 
Although in general NM can fail because of the dete- 
rioration of the simplex geometry or lack of sufficient 
decrease, the convergence of this method can be greatly 
improved by adopting a restarting strategy. Although 
other DS methods, such as the gradient descent, can per- 
form better for smooth functions, these are not robust to 
the noise which makes the objective function non-smooth 
under experimental conditions. NM has the ability to 
explore neighboring valleys with better local optima and 
likewise this exploring feature usually allows NM to over- 
come non-smoothnesses. We verified that the gradient 
descent minimization algorithm is not able to converge 
to the ground state of our Hamiltonian under the experi- 
mental conditions, mainly due to the poissonian nature of 
our photon source and the accidental counts of the detec- 
tion system, while NM converged to the global minimum 
in most optimization runs. 

Computational Scaling 

In this section, we demonstrate the polynomial scaling 
of each iteration of our algorithm with respect to system 
size, and contrast that with the exponential scaling of 
the current best-known classical algorithm for the same 
task. Suppose that the algorithm has progressed to an 
iteration j in which we have prepared a state vector l?/^-^) 
which is stored in n qubits and parameterized by the set 
of parameters {Oj}. 



We wish to find the average value of the Hamiltonian, 
{%) on this state. We will assume that there are M 
terms comprising the Hamiltonian, and assume that M 
is polynomial in the size of the physical system of interest. 
Without loss of generality, we select a single term from 
the Hamiltonian, 1-Li that acts on k bits of the state, and 
denote the average of this term by {1-Li) = h {a) where h 
is a constant and a is the /c— fold tensor product of Pauli 
operators acting on the system. As the expectation value 
of a tensor product of an arbitrary number of Pauli oper- 
ators can be measured in constant time and the spectrum 
of each of these operators is bounded, if the desired preci- 
sion on the value is given by p, we expect the cost of this 
estimation to be 0{\h\'^/p'^) repetitions of the prepara- 
tion and measurement procedure. Thus we estimate the 
cost of each function evaluation to be 0{\hmax\'^M/p^). 
For most modern classical minimization algorithms (in- 
cluding the Nelder-Mead simplex method [11 ), the cost 
of a single update step, scales linearly or at worst poly- 
nomially in the number of parameters included in the 
minimization [12]. By assumption, the number of pa- 
rameters in the set {Oj}^ is polynomial in the system 
size. Thus the total cost per iteration is roughly given 
by 0{n^\hmax\'^M/p'^) for some small constant r which 
is determined by the encoding of the quantum state and 
the classical minimization method used. 



Contrasting this to the situation where the entire al- 
gorithm is performed classically, a much different scaling 
results. Storage of the quantum state vector \iIj^) using 
currently known exact encodings of quantum states, re- 
quires knowing 2'^ complex numbers. Moreover, given 
this quantum state, the computation of the expectation 
value {a) = {i/j^ \ a \iIj^ ) using modern methods requires 
0(2^) floating point operations. Thus a single function 
evaluation is expected to require exponential resources 
in both storage and computation when performed on a 
classical computer. Moreover, the number of parameters 
which a classical minimization algorithm must manip- 
ulate to represent this state exactly is 2^. Thus per- 
forming even a single minimization step to determine 
|^j+ij> requires an exponential number of function evalu- 
ations, each of which carries an exponential cost. One 
can roughly estimate the scaling of this procedure as 
0(M2^(^+i)) 

From this coarse analysis, we conclude that our algo- 
rithm attains an exponential advantage in the cost of a 
single iteration over the best known classical algorithms, 
provided the assumptions on the Hamiltonian and quan- 
tum state are satisfied. While convergence to the final 
ground state must still respect the known complexity 
QMA-Complete complexity of this task [13], we believe 
this still demonstrates the value of our algorithm, espe- 
cially in light of the limited quantum resource require- 



ments. 

Mapping from the state parameters to the chip 
phases. The set of phases {^i}, which uniquely identi- 
fies the state \ip)^ is not equivalent to the phases which 
are written to the photonic circuit {^i}, since the chip 
phases are also used to implement the desired measure- 
ment operators cr^ (g) cr/3. Therefore, knowing the desired 
state parameters and measurement operator we compute 
the appropriate values of the chip phases on the CPU at 
each iteration of the optimization algorithm. 

EXPERIMENTAL DETAILS 

Estimation of the error on (H) 

We performed measurements of the statistical and sys- 
tematic errors that affect our computation of (H). 
Statistical errors Statistical errors due to the Poisso- 
nian nature of single photon statistics are intrinsic to the 
estimation of expectation values in quantum mechanics. 

These errors can be arbitrarily reduced at a sublinear 
cost of measurement time (i.e. efficiently) since the mag- 
nitude of error is proportional to the square root of the 
count rate. We experimentally measured the standard 
deviation of an expectation value (Hi) for a particular 
state using 50 trials. The total average coincidence rate 
was ~1500/s. The standard deviation was found to be 
37KJ/mol, which is comparable with the error observed 
in the measurement of the ground state energy shown in 
Fig. 4 

The minima of the potential energy curve was deter- 
mined by a generalized least squares procedure to fit a 
quadratic curve to the experimental data points in the re- 
gion R = (80, 100) pm, as is common in the use of trust 
region searches for minima [14 , using the inverse exper- 
imentally measured variances as weights. Covariances 
determined by the generalized least squares procedure 
were used as input to a Monte Carlo sampling procedure 
to determine the minimum energy and equilibrium bond 
distance as well as their uncertainties assuming Gaus- 
sian random error. The uncertainties reported represent 
standard deviations. Sampling error in the Monte Carlo 
procedure was 3 x 10~^ pm for the equilibrium bond dis- 
tance and 3 x 10~^ MJ/mol for the energy. 
Systematic errors In all the measurements described 
above we observed a constant and reproducible small 
shift, e = 50KJ/mol^ of the expectation value with re- 
spect to the theoretical value of the energy. There are 
at least three effects which contribute to this systematic 
error. 

Firstly, the down-conversion source that we use in our 
experiment does not produce the pure two photon state 
that is required for high-fidelity quantum interference. In 
particular, higher order photon number terms and, more 
significantly, photon distinguishability both degrade the 
performance of our entangling gate and thus the prepa- 
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ration of the state \ip). This results in a shift of the 
measured energy (?/;| H IV^). Higher order terms could be 
effectively eliminated by use of true single photon sources 
(such as quantum dots or nitrogen vacancy centers in di- 
amond), and there is no fundamental limit to the degree 
of indistinguishability which can be achieved through im- 
proved state engineering. 

Secondly, imperfections in the implementation of the 
photonic circuit also reduce the fidelity with which |?/;) 
is prepared and measured. Small deviations from de- 
signed beamsplitter refiectivities and interferometer path 
lengths, as well as imperfections in the calibration of 
voltage-controlled phases shifters used to manipulate the 
state, all contribute to this effect. However, these are 
technological limitations that can be greatly improved in 
future realizations. 

Finally, unbalanced input and output coupling effi- 
ciency also results in skewed two-photon statistics, again 
shifting the measured expectation value of {%). 

Another systematic effect that can be noted in Fig. 
4 is that the magnitude of the error on the experimen- 
tal estimation of the ground state energy increases with 
R. This is due to the fact that as R increases, the first 
and second excited eigenstates of this Hamiltonian be- 
come degenerate, resulting in increased difficulty for the 
classical minimization, generating mixtures of states that 
increases the overall variance of the estimation. 

Count rate 

In our experiment the mean count rate, which directly 
determines the statistical error, was approximately 2000- 
4000 twofold events per second. For example, for the 
bond dissociation curve we measured about 100 points 
per optimization run. The expectation value of a given 
Hamiltonian was reconstructed at each point from four 
two-qubit Pauli measurements. In the full dissociation 
curve we found the ground states of 79 Hamiltonians. 
Hence the full experiment was performed in about 158 
hours. 

State preparation is relatively fast, requiring a few mil- 
liseconds to set the phases on the chip. However 17 
seconds are required for cooling the chip, resulting in a 



duty-cycle of ^ 5%. The purpose of this is to overcome 
instability of the fibre-to-chip coupling due to thermal 
expansion of the chip during operation. This will not be 
an issue in future implementations where fibres will be 
permanently fixed to the chip's facets. Moreover the ther- 
mal phase shifters used here will also likely be replaced by 
alternative technologies based on the electro-optic effect. 
Brighter single photon sources will considerably reduce 
the measurement time. 
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